Case_InfroProp

Business Case - InfoProp

O objetivo deste trabalho é realizar uma análise minuciosa em dados imobiliários da cidade de São Paulo no ano de 2017. Através da limpeza e a vizualização dos relacionamentos entre as variáveis disponíves, a ideia é tentar captar insights importantes e, em seguida, construir um modelo de precificação baseado nas informações encontradas.

Etapa 1 - Coleta dos Dados

Para iniciar a análise, começarei pela abertura do arquivo através da função base do sistema e, em seguida, com o auxílio do pacote dplyr, tranformá-la-ei em um tbl.df para obter as melhores possibilidades de manipulação de dados.

# Carrehando pacote necessario
library(dplyr)

## 
## Attaching package: 'dplyr'

## The following objects are masked from 'package:stats':
## 
##     filter, lag

## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union

# Obtencao dos dados
imoveis <- as.tbl(read.csv2('Crawler_vivareal_venda_2017_desafio.csv', stringsAsFactors = F, encoding = "UTF-8",
                            na.strings = ""))
glimpse(imoveis)

## Observations: 7,296
## Variables: 19
## $ X.U.FEFF.date   <chr> "25/07/17", "25/07/17", "25/07/17", "25/07/17"...
## $ rua             <chr> "Coronel Artude Paula Ferreira", "Coronel Artu...
## $ numero          <chr> NA, "95", "148", "272", "450", "506", "554", "...
## $ bairro          <chr> "Vila Nova Conceição", "Vila Nova Conceição", ...
## $ condominio      <chr> NA, "Condomínio Loose In Vila Nova ", "Condomí...
## $ price           <int> 1000000, 1000000, 1000000, 1000000, 1000000, 1...
## $ area            <int> 71, 67, 45, 79, 70, 70, 48, 85, 90, 100, 73, 7...
## $ price_by_sqm    <chr> "14.084", "14.925", "22.222", "12.658", "14.28...
## $ condominium_fee <chr> "R$1.450 ", "R$1.280 ", "R$700 ", "R$648 ", "R...
## $ iptu            <chr> "R$280 ", "R$2.200 ", NA, "R$220 ", "R$240 ", ...
## $ rooms           <int> 2, 2, 1, 1, 1, 1, 1, 1, 1, 3, 2, 1, 2, 2, 3, 1...
## $ bathrooms       <int> 1, 1, 1, 1, 1, 1, 2, 2, 2, 3, 2, 2, 2, 2, 4, 1...
## $ garage          <int> 1, 1, 1, 1, 1, 2, 1, 1, 1, 1, 2, 1, 1, 1, 1, 1...
## $ agent           <chr> "Riala Empreendimentos Imobiliários", "Nalva I...
## $ agent_number    <chr> "(11) 5574-0888 ", "(11) 95772-7990 ", "(11) 9...
## $ latitude        <chr> "-235.969", "-23.596.244", "-2.359.186", "-2.3...
## $ longitude       <chr> "-466.699", "-46.669.762", "-46.679.232", "-46...
## $ property_d      <int> 74376409, 78347836, 77289482, 73325473, 728155...
## $ url             <chr> "https://www.vivareal.com.br/imovel/apartament...

Etapa 2 - Limpeza e Preparação dos dados

O primeiro passo será a limpeza dos dados. Após realizar a abertura do arquivo, identifiquei as necessidades de correção nas variáveis, das quais se destacam: problema nos encodings das variáveis data, latitude e longitude; e o formato de texto das variáveis numéricas. Os pacotes selecionados para tais tarefas foram o lubridate e stringr. Ao término do processo, percebe-se a diferença no novo dataframe.

### Data celaning
imoveis_clean = imoveis

# Pacotes necessários
library(stringr)
library(lubridate)

## 
## Attaching package: 'lubridate'

## The following object is masked from 'package:base':
## 
##     date

# data
imoveis_clean  <-  imoveis

# data
names(imoveis_clean) <- str_replace(names(imoveis_clean), 'X.U.FEFF.date', 'Data')
imoveis_clean$Data <- as.Date(imoveis$X.U.FEFF.date, format = '%d/%m/%Y')
year(imoveis_clean$Data) <- year(imoveis_clean$Data)+2000

# lat e long
imoveis_clean$longitude <- gsub(".", "", imoveis$longitude, fixed = TRUE) %>%
  sub("(\\d{2})", "\\1\\.", .) %>% 
  as.numeric()
imoveis_clean$latitude <- gsub(".", "", imoveis$latitude, fixed = TRUE) %>%
  sub("(\\d{2})", "\\1\\.", .) %>% 
  as.numeric()

# price_by_sqm
imoveis_clean$price_by_sqm <- as.numeric(imoveis_clean$price_by_sqm)

# condominium_fee e iptu
imoveis_clean$condominium_fee <- str_replace_all(imoveis$condominium_fee, 'R\\$', '') %>% 
  str_replace_all('[:punct:]', '') %>% 
  as.numeric()

imoveis_clean$iptu <- str_replace_all(imoveis$iptu, 'R\\$', '') %>% 
  str_replace_all('[:punct:]', '') %>% 
  as.numeric()

# visualizacao do novo df
glimpse(imoveis_clean)

## Observations: 7,296
## Variables: 19
## $ Data            <date> 2017-07-25, 2017-07-25, 2017-07-25, 2017-07-2...
## $ rua             <chr> "Coronel Artude Paula Ferreira", "Coronel Artu...
## $ numero          <chr> NA, "95", "148", "272", "450", "506", "554", "...
## $ bairro          <chr> "Vila Nova Conceição", "Vila Nova Conceição", ...
## $ condominio      <chr> NA, "Condomínio Loose In Vila Nova ", "Condomí...
## $ price           <int> 1000000, 1000000, 1000000, 1000000, 1000000, 1...
## $ area            <int> 71, 67, 45, 79, 70, 70, 48, 85, 90, 100, 73, 7...
## $ price_by_sqm    <dbl> 14.084, 14.925, 22.222, 12.658, 14.285, 14.285...
## $ condominium_fee <dbl> 1450, 1280, 700, 648, 1090, 1100, 1100, 800, 7...
## $ iptu            <dbl> 280, 2200, NA, 220, 240, 240, 120, 140, 100, N...
## $ rooms           <int> 2, 2, 1, 1, 1, 1, 1, 1, 1, 3, 2, 1, 2, 2, 3, 1...
## $ bathrooms       <int> 1, 1, 1, 1, 1, 1, 2, 2, 2, 3, 2, 2, 2, 2, 4, 1...
## $ garage          <int> 1, 1, 1, 1, 1, 2, 1, 1, 1, 1, 2, 1, 1, 1, 1, 1...
## $ agent           <chr> "Riala Empreendimentos Imobiliários", "Nalva I...
## $ agent_number    <chr> "(11) 5574-0888 ", "(11) 95772-7990 ", "(11) 9...
## $ latitude        <dbl> -23.59690, -23.59624, -23.59186, -23.59758, -2...
## $ longitude       <dbl> -46.66990, -46.66976, -46.67923, -46.67164, -4...
## $ property_d      <int> 74376409, 78347836, 77289482, 73325473, 728155...
## $ url             <chr> "https://www.vivareal.com.br/imovel/apartament...

summary(imoveis_clean)

##       Data                rua               numero         
##  Min.   :2017-07-25   Length:7296        Length:7296       
##  1st Qu.:2017-07-25   Class :character   Class :character  
##  Median :2017-07-25   Mode  :character   Mode  :character  
##  Mean   :2017-07-25                                        
##  3rd Qu.:2017-07-25                                        
##  Max.   :2017-07-25                                        
##                                                            
##     bairro           condominio            price         
##  Length:7296        Length:7296        Min.   :   18000  
##  Class :character   Class :character   1st Qu.: 1200000  
##  Mode  :character   Mode  :character   Median : 2900000  
##                                        Mean   : 4253352  
##                                        3rd Qu.: 5624205  
##                                        Max.   :48000000  
##                                                          
##       area          price_by_sqm     condominium_fee       iptu      
##  Min.   :   18.0   Min.   :  2.432   Min.   :     1   Min.   :    1  
##  1st Qu.:   86.0   1st Qu.: 13.571   1st Qu.:  1071   1st Qu.:  262  
##  Median :  166.0   Median : 17.026   Median :  2000   Median :  800  
##  Mean   :  215.3   Mean   : 17.709   Mean   :  2860   Mean   : 2074  
##  3rd Qu.:  293.0   3rd Qu.: 20.755   3rd Qu.:  4000   3rd Qu.: 1887  
##  Max.   :30700.0   Max.   :224.000   Max.   :460000   Max.   :50000  
##  NA's   :40        NA's   :40        NA's   :934      NA's   :1865   
##      rooms         bathrooms          garage          agent          
##  Min.   : 1.00   Min.   : 1.000   Min.   : 1.000   Length:7296       
##  1st Qu.: 2.00   1st Qu.: 2.000   1st Qu.: 1.250   Class :character  
##  Median : 3.00   Median : 3.000   Median : 3.000   Mode  :character  
##  Mean   : 2.91   Mean   : 3.451   Mean   : 2.991                     
##  3rd Qu.: 4.00   3rd Qu.: 5.000   3rd Qu.: 4.000                     
##  Max.   :10.00   Max.   :54.000   Max.   :16.000                     
##  NA's   :16      NA's   :761      NA's   :194                        
##  agent_number          latitude        longitude        property_d      
##  Length:7296        Min.   :-23.64   Min.   :-46.70   Min.   :26572520  
##  Class :character   1st Qu.:-23.60   1st Qu.:-46.67   1st Qu.:68891954  
##  Mode  :character   Median :-23.59   Median :-46.67   Median :76929697  
##                     Mean   :-23.59   Mean   :-46.67   Mean   :73681658  
##                     3rd Qu.:-23.59   3rd Qu.:-46.67   3rd Qu.:81281218  
##                     Max.   :-23.56   Max.   :-46.64   Max.   :83874638  
##                     NA's   :1661     NA's   :1661                       
##      url           
##  Length:7296       
##  Class :character  
##  Mode  :character  
##                    
##                    
##                    
## 

Etapa 3 - Análise Exploratória dos Dados

Com a preparacao dos dados concluída, pode-se prosseguir à análise exploratória. Para tal objetivo, procederemos a uma série de visualizações para ilustrar melhor a relação entre a variável de interesse em relação as demais.

## Visualizacao interativa
library(leaflet)

leaflet(data = imoveis_clean) %>% addTiles() %>% 
  addMarkers(~longitude, ~latitude, popup = ~as.character(property_d), 
             clusterOptions = markerClusterOptions()
  )

## Warning in validateCoords(lng, lat, funcName): Data contains 1661 rows with
## either missing or invalid lat/lon values and will be ignored

A começar com uma visualização interativa dos imóveis no mapa da cidade, através do pacote leaflet, consegue-se obersvar uma série de 20 coordenadas geográficas distantes da maioria das observações, configurando um ponto de interesse para se analisar com mais calma em seguida.

# visualizacao das variaveis continuas
library(ggplot2)

cont_vars <- c('area', 'iptu','condominium_fee')
plots_cont <- list()
for (var in cont_vars) {
  plots_cont[[var]] <- ggplot(imoveis_clean, aes_string(log(imoveis_clean[[var]]), log(imoveis_clean[['price']]))) + 
    geom_point() +
    geom_smooth(method = lm, se = F) + 
    ggtitle(paste(var, 'x price')) +
    theme_minimal()
  print(plots_cont[[var]])
}

## Warning: Removed 40 rows containing non-finite values (stat_smooth).

## Warning: Removed 40 rows containing missing values (geom_point).

## Warning: Removed 1865 rows containing non-finite values (stat_smooth).

## Warning: Removed 1865 rows containing missing values (geom_point).

## Warning: Removed 934 rows containing non-finite values (stat_smooth).

## Warning: Removed 934 rows containing missing values (geom_point).

Com o auxílio do pacote ggplot2, pode-se vislumbrar a relação gráfica entre as variáveis contínuas do dataframe. É interessante notar a ocorrência de outliers na relação entre ‘price’, ‘area’ e ‘condominium_fee’, gerando outras observações a serem olhadas mais de perto.

# visualizacao das variaveis factor
fact_vars <- c('bairro', 'rooms', 'bathrooms','garage')
plots_fact <- list()
for (var in fact_vars) {
  plots_fact[[var]] <- ggplot(imoveis_clean, aes_string(x = as.factor(imoveis_clean[[var]]), log(imoveis_clean[['price']]))) +
    geom_boxplot() + 
    ggtitle(paste(var, 'x price')) +
    theme_minimal()
  print(plots_fact[[var]])
}

plots_fact_c <- list()
for (var in fact_vars) {
  plots_fact_c[[var]] <- ggplot(imoveis_clean, aes_string(x = as.factor(imoveis_clean[[var]]))) + 
    geom_bar(stat = 'count') +
    ggtitle(paste('Total Observations in ',var)) +
    theme_minimal()
  print(plots_fact_c[[var]])
}

Por fim toma-se nota das relações entre as variáveis categóricas com auxílio do boxplots e gráficos de barra. Praticamente, todas as variáveis possuem algum outlier interessante, ou poucas observações de valores muito elevados, assim se fará uma sensível análise caso a caso para entender melhor suas particularidades. Outro ponto importante é a inexistência de outros bairros a não ser Vila Nova Conceição.

Etapa 4 - Outliers

Nesta seção, realizaremos um breve estudo dos outliers identificados ao longo da análise exploratória para melhor entender seus aparecimentos.

## Investigacao outliers
imoveis_clean %>% filter(log(area) > 10 | log(price) < 10 | log(condominium_fee) > 10) 

## # A tibble: 3 x 19
##   Data       rua   numero bairro condominio  price  area price_by_sqm
##   <date>     <chr> <chr>  <chr>  <chr>       <int> <int>        <dbl>
## 1 2017-07-25 <NA>  <NA>   Vila ~ <NA>       6.90e6 30700          224
## 2 2017-07-25 Coro~ <NA>   Vila ~ <NA>       1.80e4   226           79
## 3 2017-07-25 Mont~ <NA>   Vila ~ "Condomín~ 4.80e5    32           15
## # ... with 11 more variables: condominium_fee <dbl>, iptu <dbl>,
## #   rooms <int>, bathrooms <int>, garage <int>, agent <chr>,
## #   agent_number <chr>, latitude <dbl>, longitude <dbl>, property_d <int>,
## #   url <chr>

imoveis_clean %>% filter(bathrooms == 54)

## # A tibble: 1 x 19
##   Data       rua   numero bairro condominio  price  area price_by_sqm
##   <date>     <chr> <chr>  <chr>  <chr>       <int> <int>        <dbl>
## 1 2017-07-25 Balt~ <NA>   Vila ~ "Condomín~ 3.18e6   220         14.5
## # ... with 11 more variables: condominium_fee <dbl>, iptu <dbl>,
## #   rooms <int>, bathrooms <int>, garage <int>, agent <chr>,
## #   agent_number <chr>, latitude <dbl>, longitude <dbl>, property_d <int>,
## #   url <chr>

# outliers de area e price: informações muito discrepantes em relação à realidade. Provavelmente, erros de digitação
imoveis_clean$area[which(imoveis$area == 30700)] <- NA
imoveis_clean <- imoveis_clean %>% filter(!(price == 18000))
imoveis_clean$condominium_fee[which(imoveis$condominium_fee == 460000)] <- NA

# o apartamento com mais de 54 banheiros parece muito fora da media em relacao ao seu preco e area
imoveis_clean$bathrooms[which(imoveis$bathrooms == 54)] <- NA

Os valores de ‘price’, ‘area’ e ‘condominium_fee’, assim como vizualizados acima, estavam muito discrepantes da realidade sem motivo aparente, dando uma boa evidência de serem erros de digitação. Por não haver forma melhor de checar sua veracidade, e para não gerar viés indesejado na análise, optei por retirá-los da amostra. Como o único imóvel com um número muito acima da média de banheiros não seguia outros em relação a área e valor, optei pela mesma abordagem também.

Etapa 5 - Engenharia de Variáveis

Para a construção de um modelo de precificação eficiente, tem-se de extrair o máximo de informação de cada uma das variáveis e, ao mesmo tempo, evitar o overfitting, portanto a construção de novas variáveis é um passo fundamental da análise.

A começar por features já existentes, faz-se a agregação de algumas variáveis categóricas com o intuito de diminuir classes com pouca representatividade na amostra como visto nos gráficos de barra. A ideia é evitar o overfitting.

# agregação dos valores para não causar overfitting pela baixa quantidade de valores
roomsBathGarage <- imoveis_clean %>% 
  select(c(rooms, bathrooms, garage, property_d)) %>% 
  mutate(rooms_coerc = ifelse(rooms > 4, 4, rooms),
         bathrooms_coerc = ifelse(bathrooms > 6, 6, bathrooms),
         garage_coerc = ifelse(garage > 5, 5, garage)) %>% 
  select(-c(rooms, bathrooms, garage))

Como as variáveis de identificação, em seu estado bruto, não são relevantes à análise, executarei a construção de outras com possibilidade maior de fornecer informação. Após realizar sua criação, faço plots demonstrando suas informações em relação às variáveis antigas.

# criacao dos dfs comas supostas novas variaveis de interesse
imoveis_clean_rua <- imoveis_clean %>% 
  group_by(rua) %>%
  summarise(total_a = n(),
            media_p = mean(price_by_sqm, na.rm = T)) %>%
  filter(!(is.na(rua)))

imoveis_clean_condominio <- imoveis_clean %>% 
  group_by(condominio) %>%
  summarise(total_a = n(),
            media_p = mean(price_by_sqm, na.rm = T)) %>%
  filter(!(is.na(condominio))) 

# como ha agencias sem nenhuma entrada, filtro antes de achar as variaves finais
imoveis_clean_agent <- imoveis_clean %>% 
  group_by(agent) %>%
  summarise(total_a = n(),
            media_p = mean(price_by_sqm, na.rm = T)) %>%
  filter(!(is.na(agent))) %>% 
  filter(!(is.na(media_p))) 

# plot total_a
var_t <- list(imoveis_clean_rua %>% arrange(desc(total_a)), 
              imoveis_clean_condominio %>% arrange(desc(total_a)), 
              imoveis_clean_agent %>% arrange(desc(total_a)))

plots_total <- list()

for (i in 1:length(var_t)){
  plots_total[[i]] <- ggplot(var_t[[i]], aes(x = factor(var_t[[i]][[1]], levels = var_t[[i]][[1]]), 
                                             y = total_a)) + 
    geom_bar(stat="identity", width=.5) +
    coord_flip() +
    labs(title=paste('Anúncios x ', names(var_t[[i]])[[1]]), 
         x = names(names(var_t[[i]])[[1]]),
         y = 'Anúncios') 
  
  print(plots_total[[i]])
}

# plot da media_p 
var_m <- list(imoveis_clean_rua %>% arrange(desc(media_p)), 
              imoveis_clean_condominio %>% arrange(desc(media_p)), 
              imoveis_clean_agent %>% arrange(desc(media_p)))

plots_media <- list()

for (i in 1:length(var_m)){
  plots_media[[i]] <- ggplot(var_m[[i]], aes(x = factor(var_m[[i]][[1]], levels = var_m[[i]][[1]]), 
                                             y = media_p)) + 
    geom_bar(stat="identity", width=.5) +
    coord_flip() +
    labs(title=paste('Media x ', names(var_m[[i]])[[1]]), 
         x = names(names(var_m[[i]])[[1]]),
         y = 'Media') 
  
  print(plots_media[[i]])
}

# variaveis criadas
vars_rua <- imoveis_clean_rua %>% select(rua, anuncio_rua = total_a, valores_rua = media_p) 
vars_condo <- imoveis_clean_condominio %>% select(condominio, anuncio_condo = total_a, valores_condo = media_p)  
vars_agente <- imoveis_clean_agent %>% select(agent, anuncio_agente = total_a, valores_agente = media_p) 

Etapa 5 - Dataframe Final para Análise

Após terminar a construção de novas variáveis, prossegue-se às etapas de criação do dataset final para análise. Começarei pela junção das novas variáveis, e exclusão das identificadoras, além de ‘bairro’ que, como visto ao longo de toda análise exploratória, possui variância quase inexistente.

# df final para analise e com variaveis novas
imoveis_final <- imoveis_clean %>%
  left_join(roomsBathGarage) %>% 
  left_join(vars_rua) %>% 
  left_join(vars_condo) %>% 
  left_join(vars_agente) %>% 
  select(-c(Data, rua, numero, bairro, price_by_sqm, condominio, rooms, bathrooms, garage,agent, agent_number, 
            latitude, longitude, property_d, url))

## Joining, by = "property_d"

## Joining, by = "rua"

## Joining, by = "condominio"

## Joining, by = "agent"

glimpse(imoveis_final)

## Observations: 7,295
## Variables: 13
## $ price           <int> 1000000, 1000000, 1000000, 1000000, 1000000, 1...
## $ area            <int> 71, 67, 45, 79, 70, 70, 48, 85, 90, 100, 73, 7...
## $ condominium_fee <dbl> 1450, 1280, 700, 648, 1090, 1100, 1100, 800, 7...
## $ iptu            <dbl> 280, 2200, NA, 220, 240, 240, 120, 140, 100, N...
## $ rooms_coerc     <dbl> 2, 2, 1, 1, 1, 1, 1, 1, 1, 3, 2, 1, 2, 2, 3, 1...
## $ bathrooms_coerc <dbl> 1, 1, 1, 1, 1, 1, 2, 2, 2, 3, 2, 2, 2, 2, 4, 1...
## $ garage_coerc    <dbl> 1, 1, 1, 1, 1, 2, 1, 1, 1, 1, 2, 1, 1, 1, 1, 1...
## $ anuncio_rua     <int> 320, 320, 49, 570, 2, 94, 383, 276, 276, 383, ...
## $ valores_rua     <dbl> 18.31858, 18.31858, 19.19827, 18.79077, 15.164...
## $ anuncio_condo   <int> NA, 29, 7, 87, NA, 21, 85, 53, 53, 10, 18, NA,...
## $ valores_condo   <dbl> NA, 15.17407, 22.18257, 19.48547, NA, 16.90986...
## $ anuncio_agente  <int> 1, 3, 2, 72, 9, 73, 49, 49, 2, 19, 3, 5, 74, 1...
## $ valores_agente  <dbl> 14.08400, 13.78067, 21.11100, 18.12724, 16.408...

Em seguida, prosseguir-se-á “fatorização” do dataframe, ou seja, transformação das variáveis categóricas para o formato factor

# Função mais eficiente para featurização
imoveis_fact <- imoveis_final %>%
  mutate_at(
    .vars = vars('rooms_coerc', 'bathrooms_coerc','garage_coerc'),
    .funs = funs(as.factor(.))
  )

glimpse(imoveis_fact)

## Observations: 7,295
## Variables: 13
## $ price           <int> 1000000, 1000000, 1000000, 1000000, 1000000, 1...
## $ area            <int> 71, 67, 45, 79, 70, 70, 48, 85, 90, 100, 73, 7...
## $ condominium_fee <dbl> 1450, 1280, 700, 648, 1090, 1100, 1100, 800, 7...
## $ iptu            <dbl> 280, 2200, NA, 220, 240, 240, 120, 140, 100, N...
## $ rooms_coerc     <fct> 2, 2, 1, 1, 1, 1, 1, 1, 1, 3, 2, 1, 2, 2, 3, 1...
## $ bathrooms_coerc <fct> 1, 1, 1, 1, 1, 1, 2, 2, 2, 3, 2, 2, 2, 2, 4, 1...
## $ garage_coerc    <fct> 1, 1, 1, 1, 1, 2, 1, 1, 1, 1, 2, 1, 1, 1, 1, 1...
## $ anuncio_rua     <int> 320, 320, 49, 570, 2, 94, 383, 276, 276, 383, ...
## $ valores_rua     <dbl> 18.31858, 18.31858, 19.19827, 18.79077, 15.164...
## $ anuncio_condo   <int> NA, 29, 7, 87, NA, 21, 85, 53, 53, 10, 18, NA,...
## $ valores_condo   <dbl> NA, 15.17407, 22.18257, 19.48547, NA, 16.90986...
## $ anuncio_agente  <int> 1, 3, 2, 72, 9, 73, 49, 49, 2, 19, 3, 5, 74, 1...
## $ valores_agente  <dbl> 14.08400, 13.78067, 21.11100, 18.12724, 16.408...

Etapa 6 - Modelos Preditivos

Chega-se à parte mais relevante do trabalho: construção do modelo preditivo para precificação de imóveis. O pacote escolhido para a tarefa foi o “mlr”. Atualmente, um dos pacotes mais utilizados por usuários da linguagem, possui muitas opções de customização e facilidade de aplicação.

A primeira etapa consistirá na preparação dos dados através da padronização das variáveis contínuas, o tratamento dos valores missing e o processo conhecido como “One Hot Encoding”, ou seja, transformação das variáveis factor em dummies por levels.

## Pacote MLR

# Pacote
library(mlr)

## Loading required package: ParamHelpers

set.seed(666)

# Imputador dos valores missing
imp <- impute(
  imoveis_fact,
  classes = list(
    factor = imputeMode(),
    integer = imputeMean(),
    numeric = imputeMean()
  )
)
imoveis_imp <- imp$data

# Sumarizando as colunas
summarizeColumns(imoveis_imp) %>%
  knitr::kable(digits = 2) # gerador de tabelas muito prático do knitr

name	type	na	mean	disp	median	mad	min	max	nlevs
price	numeric	0	4253932.96	4780029.51	2900000.00	2824353.00	250000.00	4.800e+07	0
area	numeric	0	211.11	158.71	168.00	133.43	18.00	1.683e+03	0
condominium_fee	numeric	0	2860.32	5790.93	2450.00	1868.08	1.00	4.600e+05	0
iptu	numeric	0	2073.65	3929.10	1500.00	1156.43	1.00	5.000e+04	0
rooms_coerc	factor	0	NA	0.63	NA	NA	1108.00	2.698e+03	4
bathrooms_coerc	factor	0	NA	0.72	NA	NA	899.00	2.034e+03	6
garage_coerc	factor	0	NA	0.73	NA	NA	1221.00	1.970e+03	5
anuncio_rua	numeric	0	245.57	137.24	245.57	113.52	1.00	5.700e+02	0
valores_rua	numeric	0	17.65	3.34	17.65	1.75	6.25	4.866e+01	0
anuncio_condo	numeric	0	38.53	22.09	38.53	0.00	1.00	1.180e+02	0
valores_condo	numeric	0	17.49	3.95	17.49	0.00	4.20	5.177e+01	0
anuncio_agente	numeric	0	42.37	33.57	38.00	37.06	1.00	1.350e+02	0
valores_agente	numeric	0	17.70	2.81	18.02	2.10	2.43	5.839e+01	0

# Normalização
imoveis_norm <- normalizeFeatures(imoveis_imp, target = "price")

# One hot encoding
imoveis_preProc <- createDummyFeatures(
  imoveis_norm, target = "price",
  cols = c(
    "rooms_coerc",
    "bathrooms_coerc",
    "garage_coerc"
  )
)

Com os dados preparados para análise, cria-se o objeto básico utilizado neste pacote: task (no caso, uma task específica para regressão). A partir dele, faz-se a separação hold-out (treino x teste)

# Criação do Task para Regressão (objeto usado nas operações do pacote)
task <- makeRegrTask(id = 'Imoveis_SP', data = imoveis_preProc, target = 'price')

# train test split (holdout é o método de resample de separar em train/test)
holdout <- makeResampleInstance('Holdout', task)
task_train <- subsetTask(task, holdout$train.inds[[1]])
task_test <- subsetTask(task, holdout$test.inds[[1]])

Com os objetos criados, pode-se iniciar o treinamento dos modelos. Serão usados 2 neste trabalho: um linear simples sem qualquer tipo de tratamento, para se usar como benchmark, e um XGBoost usando validações cruzadas e tunagem de hiper-parâmetros. A escolha do modelo deveu-se pelas relações não lineares entre as variáveis númericas, e seu sistema de uso de gradiente para minimização dos resíduos. 5 validações cruzadas serão utilizadas para evitar o overfitting e a série de hiper-parâmetros escolhidos para tunagem podem ter sua explicação melhor desenvolvida em https://github.com/dmlc/xgboost.

## Modelo linear básico

# modelo padrão linear
regr.lm <- makeLearner(id = 'lm', 'regr.lm')

# treinamento modelo
lm <- train(regr.lm, task_train)

## XGBoost

# Criação do modelo simples
xgb_model <- makeLearner("regr.xgboost")

## Warning in makeParam(id = id, type = "numeric", learner.param = TRUE, lower = lower, : NA used as a default value for learner parameter missing.
## ParamHelpers uses NA as a special value for dependent parameters.

# Grid de alguns parâmetros
xgb_params <- makeParamSet(
  # The number of trees in the model (each one built sequentially)
  makeIntegerParam("nrounds", lower = 100, upper = 2000),
  # number of splits in each tree
  makeIntegerParam("max_depth", lower = 1, upper = 10),
  # "shrinkage" - prevents overfitting
  makeNumericParam("eta", lower = .1, upper = .5),
  # L2 regularization - prevents overfitting
  makeNumericParam("lambda", lower = -1, upper = 0, trafo = function(x) 10^x)
)

# Controle da tunagem randômica
control <- makeTuneControlRandom(maxit = 1) # maxit represente tempo máximo de 1 min para iterações

# Plano de resampling
resample_desc <- makeResampleDesc("CV", iters = 5)

# tunador dos hiper parâmetros
tuned_params <- tuneParams(
  learner = xgb_model,
  task = task_train,
  resampling = resample_desc,
  measures = rmse,       # R-Squared performance measure, this can be changed to one or many
  par.set = xgb_params,
  control = control
)

## [Tune] Started tuning learner regr.xgboost for parameter set:

##              Type len Def       Constr Req Tunable Trafo
## nrounds   integer   -   - 100 to 2e+03   -    TRUE     -
## max_depth integer   -   -      1 to 10   -    TRUE     -
## eta       numeric   -   -   0.1 to 0.5   -    TRUE     -
## lambda    numeric   -   -      -1 to 0   -    TRUE     Y

## With control class: TuneControlRandom

## Imputation value: Inf

## [Tune-x] 1: nrounds=1774; max_depth=3; eta=0.465; lambda=0.374

## [Tune-y] 1: rmse.test.rmse=1209875.5930818; time: 0.3 min

## [Tune] Result: nrounds=1774; max_depth=3; eta=0.465; lambda=0.374 : rmse.test.rmse=1209875.5930818

# Modelo tunado
xgb_tuned_model <- setHyperPars(
  learner = xgb_model,
  par.vals = tuned_params$x
)

# Verificação da performance do modelo tunado usando as validações cruzadas usadas para tunagem
resample(xgb_tuned_model,task_train,resample_desc,measures = list(rmse, rsq))

## Resampling: cross-validation

## Measures:             rmse      rsq

## [Resample] iter 1:    1009558.34928650.9477924

## [Resample] iter 2:    1616862.04355040.8780969

## [Resample] iter 3:    1167638.40041130.9387638

## [Resample] iter 4:    1076768.89157160.9451954

## [Resample] iter 5:    1423415.55945190.9258245

## 

## Aggregated Result: rmse.test.rmse=1279247.7057512,rsq.test.mean=0.9271346

## 

## Resample Result
## Task: Imoveis_SP
## Learner: regr.xgboost
## Aggr perf: rmse.test.rmse=1279247.7057512,rsq.test.mean=0.9271346
## Runtime: 20.0695

# Treinamento
XGBoost <- train(xgb_tuned_model, task_train)

Etapa 7 - Visualização dos Resultados

Agora, com os modelos devidamente treinadss, podemos fazer as comparações dos resultados, quais variáveis possuem maior efeito e como é a sua influência nas previsões. A começar pela performance nos dados de teste, faz-se uma comparação de 2 indicadores de qualidade das previsões juntamente a uma ilustrativa.

# predição nos dados de teste
pred_lm <- predict(lm, task_test)

## Warning in predict.lm(.model$learner.model, newdata = .newdata, se.fit =
## FALSE, : prediction from a rank-deficient fit may be misleading

performance(pred_lm, measures = list(rmse, rsq))

##         rmse          rsq 
## 3.124175e+06 5.905637e-01

pred_xgb <- predict(XGBoost, task_test)
performance(pred_xgb, measures = list(rmse, rsq))

##         rmse          rsq 
## 1.196479e+06 9.399483e-01

## Análise dos dados fitted e feature importance
comparacao <- data_frame(pred_lm = pred_lm[['data']][['response']], pred_xgb = pred_xgb[['data']][['response']],
                         price = getTaskData(task_test)[['price']])

# Plot predictions (on x axis) vs actual bike rental count
comparacao %>%
  tidyr::gather(tipo, valor, price, pred_lm, pred_xgb) %>%
  filter(!(tipo == 'pred_xgb')) %>% 
  ggplot(aes(valor, fill = tipo)) +
  ggtitle('LM') +
  geom_density(alpha = .5)

comparacao %>%
  tidyr::gather(tipo, valor, price, pred_lm, pred_xgb) %>%
  filter(!(tipo == 'pred_lm')) %>% 
  ggplot(aes(valor, fill = tipo)) +
  ggtitle('LM') +
  geom_density(alpha = .5)

Como demonstrado acima, e já esperado, o modelo xgbost apresentou índices superiores nos indicadores, além de, visualmente, ter uma aderência maior aos dados de teste. Portnato, para encerramento da análise, usarei o pacote ‘iml’, recomendado pelos criadores do ‘mlr’, para aumentar a compreensão desse modelo mais complexo. Primeiro, descobrindo qual a importância de cada uma das features e, em seguida, como os valores da mais importante influenciam os valores previstos pelo modelo.

# Pacote
library(iml)
# usando Predictor$new() para criar um plot das features mais importantes
X = getTaskData(task_train)[getTaskFeatureNames(task_train)]
predictor = Predictor$new(XGBoost, data = X, y = getTaskData(task_train)['price'])
imp = FeatureImp$new(predictor, loss = "rmse")
plot(imp)

# Plot de como a variável mais importante afeta as previsões na média
pdp_area = Partial$new(predictor, feature = "area")

## Warning: The FeatureEffect class replaces the Partial class. Partial will
## be removed in future versions.

plot(pdp_area)

FIM